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ABSTRACT 

A classic method for computing the mass function of dark matter halos is provided by excursion set the- 
ory, where density perturbations evolve stochastically with the smoothing scale, and the problem of computing 
the probability of halo formation is mapped into the so-called first-passage time problem in the presence of 
a barrier. While the full dynamical complexity of halo formation can only be revealed through A^-body sim- 
ulations, excursion set theory provides a simple analytic framework for understanding various aspects of this 
complex process. In this series of paper we propose improvements of both technical and conceptual aspects of 
excursion set theory, and we explore up to which point the method can reproduce quantitatively the data from 
A^-body simulations. In paper I of the series we show how to derive excursion set theory from a path integral 
formulation. This allows us both to derive rigorously the absorbing barrier boundary condition, that in the usual 
formulation is just postulated, and to deal analytically with the non-markovian nature of the random walk. Such 
a non-markovian dynamics inevitably enters when either the density is smoothed with filters such as the top-hat 
filter in coordinate space (which is the only filter associated to a well defined halo mass) or when one considers 
non-Gaussian fluctuations. In these cases, beside "markovian" terms, we find "memory" terms that reflect the 
non-markovianity of the evolution with the smoothing scale. We develop a general formalism for evaluating 
perturbatively these non-markovian corrections, and in this paper we perform explicitly the computation of the 
halo mass function for gaussian fluctuations, to first order in the non-markovian corrections due to the use of a 
tophat filter in coordinate space. 

In paper II of this series we propose to extend excursion set theory by treating the critical threshold for 
collapse as a stochastic variable, which better captures some of the dynamical complexity of the halo formation 
phenomenon, while in paper III we use the formalism developed in the present paper to compute the effect of 
non-Gaussianities on the halo mass function. 

Subject headings: cosmology:theory — dark matter:halos — large scale structure of the universe 



1. INTRODUCTION 

The computation of the mass function of dark matter ha- 
los is a central problem in modern cosmology. In par- 
ticular, the high-mass tail of the distribution is a sensitive 
probe of primordial non-Gaussianities ([Matarrese et al. 1986; 
Moscardini et al. 1991; Kovama et al. 1999; Matarres e et al.l 
M)0; Robi nson & Baker 2000; Robinson et al. 2000). Var- 
ious planned large-scale galaxy surveys, both ground based 
(DES, PanSTARRS and LSST) and on satellite (EUCLID and 
ADEPT) can detect the effect of primordial non-Gaussianities 
on the mass d istributio n of dark ma tter halos (see e.g. 
iDalal et a l.1 (12008); Carbone et al.l 12008)). Of course, this 
also requires reliable theoretical predictions for the mass func- 
tion, first of all when the primordial fluctuations are taken 
to be gaussian, and then including non-Gaussian corrections. 
Furthermore, the halo mass function is both a sensitive probe 
of cosmological parameters and a crucial ingredient when one 
studies the dark matter distribution, as well as the formation, 
evolution and distribution of galaxies, so its accurate predic- 
tion is obviously important. 

The formation and evolution of dark matter halos is a highly 
complex dynamical process, and a detailed understanding of it 
can only come through large-scale A^-body simulations. Some 
analytical understanding is however also desirable, both for 
obtaining a better physical intuition, and for the flexibility un- 



der changes of models or parameters (such as cosmological 
model, shape of the non-Gaussianities, etc.) that is the advan- 
tage of analytical results over very timing consuming numer- 
ical simulations. 

Analytic techniques generally start by modelling the col- 
lapse as spherical or ellipsoidal. However, A^-body simula- 
tions show that the actual process of halo formation is not 
ellipsoidal, and in fact is not even a collapse, but rather 
a messy mixture o f violent encounters , smooth accretion 
and fragmentation (ISpringel et al.l l2005h . In spite of this, 
analytical te chniques based on Press-Schecther (PS) theory 
(iPress & Sc hechter 1974) and its extension known as excur- 
sion set theory (Peac ock & Hea vens 19901: iBond et anil991h 
are able to reproduce, at least qualitatively, several proper- 
ties of dark matter halos such as their conditional and un- 
conditional mass functi on, halo accretion histories, merger 
rates and halo bias (see IZentned (l2007h for a recent review). 
However, at the quantitative level, already for gaussian fluc- 
tuations the prediction of excursion set theory for the mass 
function deviate significantly from the results of A^-body sim- 
ula tions. The halo m ass function dn/dM can be written 
as dJenkins et al.ll2001h 



dM 



M2 c/lnM 



(1) 
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where n(M) is the number density of dark matter halos of mass 
Universite de Geneve, 24 quai M, is the variance of the linear density field smoothed on a 

scale R corresponding to a mass M, and p is the average den- 
sity of the universe. In excursion set theory within a spherical 
collapse model the function /(ct) is predicted to be 



2 




0.30 



0.50 



0.70 



1.00 



1.50 2.00 



-1 3.00 



Fig . 1 . — A log-log plot of the function f{cr). The (blue) solid curve is the 
PS prediction /ps(a"). The three almost indistinguishable dashed lines are the 
Sheth-Tormen fit to the GIF simulation of Kauffmann et al. ( 1999), and the fit 
to the W-body simulations of iPillepich et alj ([2008) and Warren et al.j2003)- 
The fitting functions are summarized in Table 3 of Pillepich et al. ( 2008ft . 



./ps(a)=(- ^, 



(2) 



where Sc — 1.686 is the critical value for collapse in the 
spherical collapse model. This result can be extended to 
arbitrary redshift z reabsorbing the evolution of the vari- 
ance into Sc, so that 6c in the above result is replaced 
by Sciz) = Sc{0)/D{z), where D{z) is the linear growth fac- 
tor. This prediction can be com pared with th e exist - 

ing A^-body simulations (see e.g; Jenkins et al.l ( 2001 ): 

Warren et al. (2006); Lukic et all (120071) : iTinker et al.1 hoO^ : 
Pillepich et al.. (.2008.) : .Robertson et al. and references 

therein). The results of these simulations have been repre- 
sented by various fitting functions, see e.g. ' Sheth & TormenI 
(|l999), Sheth et al. (2001). In Fig.[T]we compare the function 
/ps((t) given in eq. to various fits to A^-body simulations, 
plotting the result against cr"'. High masses correspond to 
large smoothing radius R, i.e. low values of a and large cr"', 
so mass increases from left to right on the horizontal axis. One 
sees that the A^-body simulations are quite consistent among 
them, and that PS theory predicts too many low-mass halos, 
roughly by a factor of two, and too few high-mass halos: at 
cr~' = 3, PS theory is already off by a factor 0(10). The pri- 
mordial non-Gaussianities can be constrained by probing the 
statistics of rare events, such as the formation of the most mas- 
sive objects, so it is particularly important to model accurately 
the high-mass part of the halo mass function, first of all at the 
gaussian level. It makes little sense to develop an analytic 
theory of the non-Gaussianities, by perturbing over a gaus- 
sian theory that in the interesting mass range is already off by 
one order of magnitude. 

When searching for the origin of this failure of excursion set 
theory, one can divide the possible concerns into two classes: 

(i) Even if one accepts as a physical model for halo for- 
mation a spherical (or ellipsoidal) collapse model, there are 
formal mathematical problems in the implementation of ex- 
cursion set theory that leads to eq. 

(ii) The physical model itself is inadequate, since a spher- 
ical or even elliptical collapse model is an oversimplification 
of the actual complex process of halo formation. 

Concerning point (i), it is well known that the original ar- 
gument of Press and Schechter miscounts the number of viri- 
alized objects because of the so-called "cloud-in-cloud" prob- 



lem. In the spherical collapse model one assumes that a region 
of radius R, with a smoothed density contrast 6(R), collapses 
and virializes once 6{R) exceeds a critical value 6c — 1.686."^ 
Within PS theory, for gaussian fluctuations the distribution 
probabiUty for the density contrast is 



nps((5,5): 



1 



-57(25) 



where 



S(R) = aHR)={6Hx,R))., 



(3) 



(4) 



and the fractional volume of space occupied by virialized ob- 
jects larger than R is identified with 



Fps(R) ■■ 



1 

£/(5nps((5,5(/;))= -erfc 



V2 



(5) 



where viR) = 6c/<j{R)- As remarked already by Press and 
Schechter, this expression cannot however be fully correct. 
In fact, in the hierarchical models that we are considering the 
variance S{R) diverges as — > 0, so all the mass in the uni- 
verse must finally be contained in virialized objects. Thus, 
we should have Fps(O) = 1, while eq. (|5]) gives Fps(O) = 1/2. 
Press and Schechter corrected this simply adding by hand an 
overall factor of two. 

The reason for this failure is that the above procedure 
misses the cases in which, on a given smoothing scale R, 5{R) 
is below the threshold, but still it happened to be above the 
threshold at some scale R' > R. Such a configuration corre- 
sponds to a virialized object of mass M' > M. However, it is 
not counted in Fps{R) since on the scale S it is below thresh- 
old. Thus eq. (|5T l cannot be fully correct. 

In 'Bond et alJ (11991 ) this problem was solved by mapping 
the evolution of 6 with the smoothing scale into a stochastic 
problem. Using a sharp ^-space filter, they were able to for- 
mulate the problem in terms of a Langevin equation with a 
Dirac-delta noise. In other words, the smoothed density per- 
turbation 6 suffers a markovian stochastic motion under the 
influence of a gaussian white noise, with the variance S= cr^ 
playing the role of a time variable. In this formulation, the 
halo is defined to be formed when the smoothed density per- 
turbation 6 reaches the critical value Sc for the first time. The 
problem is therefore reduced to a "first-passage problem", 
which i s a classical subject in the theory of stochastic pro- 
cesses (Redner 2001). One may write a Fokker-Planck equa- 
tion describing the probability n((5,5') that the density pertur- 
bation acquires a given value S at a given "time" S, supple- 
mented by the absorbing barrier boundary condition that the 
probability vanishes when 6 = Sc- The solution reproduces 
eq. (|2]i, including the factor of two that Press and Schechter 
were forced to introduce by hand.^ 

However, this procedure still raises some technical ques- 
tions, that will be reviewed in more detail in Section |2] In 
short, there are two issues that deserve a deeper scrutiny. 
First, the "absorbing barrier" boundary condition n((5f ,5) = 
is a natural one, but still it is something that is imposed by 
hand, and in this sense it is really an ansatz. In the litera- 
ture for stochastic processes it is well-known that, in general, 

* More precisely, 5c has a slight dependence on the cosmological model, 
and 5c = 1 .686 is the value for a fi^ = 1 cosmology (Lacev & Cole 1993). For 
a model with + H a = 1 this dependence is computed in , Eke et al.i t l99d) . 
For Qm — 0.3, 5c is between 1.67 and 1.68, see their Fig. 1. This difference 
is however much smaller than other uncertainties in our computation. 

^ The work of Epstein ( 1983) also solves the cloud in cloud problem and 
recovers the coiTect factor of two, though the process considered therein uses 
Poisson seeds for structure formation. 
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the probability does not satisfy any simple boundar y condi- 
tion (Ivan Kampen & Oppenheimlll972l:lK"nessllll986l) . This is 
due to the fact that, when one works with a discretized time 
step, a stochastic trajectory can exit a given domain by jump- 
ing over the boundary without hitting it, unlike a continuous 
diffusion process which has to hit the boundary to exit the 
domain. Particular care must therefore be devoted to the pas- 
sage from the discrete to the continuum. As we will see, the 
passage from a discrete to a continuum formulation is indeed 
highly non-trivial when a generic filter and/or non-Gaussian 
perturbations are used. 

A second related concern is that the derivation of Bond et 
al. only works for a sharp ^-space filter. However, as we 
review in Section |2l there is no unambiguous way of asso- 
ciating a mass to a region of size R smoothed with a sharp 
^-space filter. The only unambiguous way of associating a 
mass M to a smoothing scale R is using a sharp filter in x- 
space, proportional to (R-r), in which case one has the obvi- 
ous relation M = (4/3)ttR^p. This is also the relation used in 
numerical simulations. As soon as one uses a different filter 
(such as the tophat in real space), the Langevin equation with 
gaussian Dirac-delta noise, that describes a simple markovian 
process, is replaced by a very complicated non-markovian 
dynamics dictated by a colored noise. The system acquires 
memory properties and the probability n((5,5) no longer sat- 
isfy a simple diffusion equation such as the Fokker-Planck 
equation. The same is true if the density perturbation is non- 
Gaussian. Furthermore, the correctness of the "absorbing bar- 
rier" boundary condition is now far from obvious. These dif- 
ficulties are well-known in the statistical physics community, 
where progress in solving the first-passage problem in the 
presence of a non-rnarkovian dynamics has been very limite d 
(iHanggi etalJ [19811 iWeiss etalJfT98l Ivan KampenI [T998h . 
From these considerations one concludes that the rather com- 
mon procedure of taking the analytical results of iBond et al.l 
11991), valid for a sharp filter in momentum space, and apply- 
ing them to generic filters is incorrect.^ 

These issues become even more important when one con- 
siders the evolution with smoothing scale of non-gaussian 
fluctuations, since non-gaussianities induce again a non- 
markovian dynamics, and furthermore it is important to disen- 
tangle the physically interesting non-markovian contribution 
to the halo mass function due to primordial non-gaussianities, 
from the non-markovian contribution due to the filter func- 
tion. 

Concerning point (ii) above, it is important to stress once 
again that excursion set theory is just a simple mathemat- 
ical model for a complex dynamical process. Treating the 
collapse as ellipsoi dal rather than spherical gives a rnore re - 
aHstic description (iSheth & Tormenlll999l:ISheth et al.ll2001l) . 
However, as we already mentioned, dark matter halos grow 
through a mixture of smooth accretion, violent encounters and 
fragmentations, and modeling halo collapse as spherical, or 
even ellipsoidal, is certainly an oversimplification. In addi- 
tion, the very definition of what is a dark matter halo, both in 
A^-body simulations and obser vationally, is a diffic ult problem 
(for cluster observations, see Jeltema et al.l (|2005) and refer- 
ences therein), that we will discuss in more detail in paper II. 

In this series of paper we examine systematically the above 

^ Similarly, even if the mathematical problem of solving the Fokker- 
Planck equation with a moving barrier is amenable to an elegant formulation 
{Zhang & Hui 2006), its application to the halo mass function suffers from 
the problem that for a general filter it is incorrect to assume that the probabil- 
ity n((5,S) evolves according to the Fokker-Planck equation. 



issues. In the present paper we start from excursion set the- 
ory in its simpler physical implementation, i.e. coupled to 
a spherical collapse model, and within this framework we put 
the formalism on firmer mathematical grounds. We show how 
to formulate the mathematical problem exactly in terms of a 
path integral with boundaries and particular care will be de- 
voted to the passage from the discrete to the continuum. This 
formalism allows us to obtain a number of result: first, when 
we restrict to gaussian fluctuations and sharp A:-space filter, 
in the continuum limit we recover the usual formulation of 
excursion set theory, but in this case the absorbing barrier 
boundary condition emerges automatically from the formal- 
ism, without the need of imposing it by hand. For different 
filters the problem becomes much more complicated, and we 
have to deal with a non-markovian dynamics. We will see 
that, for a generic filter, the zeroth-order term in an expan- 
sion of the non-markovian contributions gives back eq. Q, 
where cr^ is now the variance computed with the generic fil- 
ter. We then show how the non-markovian contributions can 
be computed perturbatively using our path integral formula- 
tion, and we compute explicitly, to first perturbative order, 
the halo mass function for a tophat filter in coordinate space. 
We find that the non-markovian contributions do not allevi- 
ate the discrepancy with A^-body simulations. On the con- 
trary, in the relevant mass range the full halo mass function 
is everywhere slightly lower than the one obtained from the 
markovian contribution, so in the large mass regime this cor- 
rection goes in the wrong direction. This result will not be 
a surprise to the expert reader. Already in their classical pa- 
per. Bond et al. computed the result with a tophat filter in 
coordinate space using a Monte Carlo (MC) realization of the 
trajectories obtained from a Langevin equation with colored 
noise, and found indeed that one has fewer high mass objects. 
More recently, a MC simulation of this kind has been done in 
iRobertson et al.l ( 12008 ). and our analytical result to first order 
is in agreement with their findings. 

In paper II of this series, motivated by the physical limita- 
tions of the spherical or ellipsoidal collapse model, we pro- 
pose that some of the physical complications of the realistic 
process of halo formation and growth can be included in the 
excursion set framework, at least at an effective level, by as- 
suming that the critical value for collapse is not a fixed con- 
stant Sc, as in the spherical collapse model, nor a fixed func- 
tion of the variance a^, as in the ellipsoidal collapse model, 
but rather is itself a stochastic variable, whose scattering re- 
flects a number of complicated aspects of the underlying dy- 
namics. 

Finally, in paper III of this series we apply the formalism 
developed in the present paper, together with the diffusing 
barrier model developed in paper II, to the computation of 
the halo mass function in the presence of non-Gaussian fluc- 
tuations. 

This paper is organized as follows. In Section [2 ] we rev iew 
the excursion set theory developed in iBond et al.l (119911) : in 
Section[3]we present the path integral approach to a stochastic 
problem in the presence of a barrier. In Sections|4]we special- 
ize to the cases of a sharp filter in momentum space, while in 
Section |5] we consider a generic filter In particular, in Sec- 
tion |5] we show how to deal with the non-markovian correc- 
tions to the halo mass function. Some technicalities regarding 
the delicate passage from the discrete to the continuum are 
contained in Appendices A and B. 
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2. THE COMPUTATION OF THE HALO MASS FUNCTION AS 
A STOCHASTIC PROBLEM 

The computation of the halo mass function can be formu- 
lated in terms of a st ochastic process, as is well known since 
the classical work of lBond et alJ (1199 ll) . Let us recall the pro- 
cedure, in order to set the notation and to highlight some deli- 
cate points, in particular related to the choice of the filter func- 
tion, that are important in the following. The expert reader 
might wish to move directly to Section[3] 

One considers the density contrast S{x) = [p{x) - p]/ p, 
where p is the mean mass density of the universe and x is the 
comoving position, and smooths it on some scale R, defining 

6ix,R)= J d^x'Wi\x-x'\,R)6ix'), (6) 

with a filter function W(\x-x'\,R). We denote by W(li,R) its 
Fourier transform. A simple choice is a sharp filter in A:-space, 

Wsh.rp-kik,kf) = d{kf-k), (7) 

where kf = l/R, k = |k| and 6 is the step function. Other 
common choices are a sharp filter in x-space, Wsharp-x('",^) = 

[3/(4-nR^)]e(R-r), or a gaussian filter, %^u(k,R) = e""'*'/^ 
Writing eq. (|6]l in terms of the Fourier transform we have 



S(x,R) = 



■Sik)W(k,R)e 



-ikx 



(8) 



where k= |k| . We focus on the evolution of (5(x, R) with Rata 
fixed value of x, that we can choose without loss of generality 
as X = 0, and we write S{x = Q,R) simply as 5(R). Taking the 
derivative of eq. dHJ with respect to R we get 



dS(R) 



where 



OR 

d^k ~^^,dW(k,R) 



■S(k)- 



dR 



(9) 



(10) 



Since the modes ^(k) are stochastic variables, (^(R) is a 
stochastic variable too, and eq. ^ has the form of a Langevin 
equation, with R playing the role of time, and ((R) playing 
the role of noise. When 6(R) is a gaussian variable, only its 
two-point connected correlator is non-vanishing. In this case, 
we see from eq. ( fTOl i that also ( is gaussian. The two-point 
function of 6 defines the power spectrum P(k), 



(5(k)^(k')) = (27r)^(5o(k + k')P(fc) . 
From this it follows that 



(11) 



d(\nk) A\k) 



dW(k,Ri) dW(k,R2) 



dRi 



dRi 



(12) 

where, as usual, A^(k) = k^P(k)/(2n^). For a generic filter 
function the right-hand side is a function of Ri and R2, differ- 
ent from a Dirac delta SoiRi-Ri)- In the literature on stochas- 
tic processes this case is known as colored gaussian noise. 
Things simplify considerably for a sharp fc-space filter. Us- 
ing kf = 1 /R instead of R, and defining Q{kf) = -(1 /^f)C(A;/?), 
eqs. |9|l and ( fT2] i become 



dSikf) 
dlnkf 



= Q(kF), 



and 



{Q(kFl)Q(kF2)) = AHkFl)SD(lnkFi-lnkF2) . 



(13) 



(14) 



Therefore, we have a Dirac delta noise. We can write these 
equations in an even simpler form using as "pseudotime" vari- 
able the variance S defined in eq. Using eq. (O 



(15) 



S(R)= / d(\nk)A\k)\W(k,R) 

J— 00 

For a sharp k-spacQ filter, S becomes 

n\nkf 

S{kF)= / d{lnk)A^(k), 



so 



dS 



= A\kf). 



d\nkf 

Thus, redefining finally rjikF) = Qikp)/ A^{kF), we get 



d5(S) 
dS 



= V(S), 



with 



{■n(Si)r](S2)) = S(S,-S2). 



(16) 
(17) 

(18) 
(19) 

which is a the Langevin equation with Dirac-delta noise, with 
S playing the role of time. In hierarchical power spectra, at 
R = oowe have 5 = 0, and S increases monotonically as R de- 
creases. Therefore we can start from R = 00, corresponding 
to "time" 5 = 0, where 6 = 0, and follow the evolution of S{S) 
as we decrease R, i.e. as we increase S. The fact that this 
evolution is governed by the Langevin equation means that 
S(S) performs a ra ndom walk , with respect to the "time" vari- 
able S. Following lBond et alj (Il99ll) . we refer to the evolution 
of (5 as a function of 5 as a "trajectory". In the spherical col- 
lapse model, a virialized object forms as soon as the trajectory 
exceeds the threshold 6 = 5c. In this language, the "cloud-in- 
cloud" problem of PS theory is associated with trajectories 
that make multiple crossings of the threshold, such as that 
shown in Fig.|2] If we compute the probability distribution at 
S = S2 as in PS theory, i.e. using eq. (|5]l, this trajectory does 
not contribute to Fps{R) since at this value of S it is below 
threshold. However, it has already gone above threshold at 
an earlier time Si, corresponding to a radius Ri, so it gives a 
virialized object of mass M{R\) > M(R2). This virialized ob- 
ject has been lost in Fps(^2) evaluated through eq. (|5]), in spite 
of the fact that this formula was supposed to count all objects 
with mass greater then M(R2). 

To cure the "cloud-in-cloud" problem we must consider the 
lowest value of S (or, equivalently, the highest value of R) 
for which the trajectory pierces the threshold. Similar prob- 
lems are known in statistical physics as "first-passage time" 
problems. After that, a virialized object forms and this trajec- 
tory should be excluded from further consideration. We there- 
fore consider an ensemble of trajectories, all starting from 
the initial value (5 = at initial "time" 5 = 0, and we com- 
pute the function H{5, S) that gives the probability distribu- 
tion of reaching a value 5 at "time" S. As is well known, if 
a stochastic process obeys the Langevin equation ( fTSb with a 
Dirac delta noise ( fT9l ). the corresponding distribution function 
is a solution of the Fokker-Planck (FP) equation, 

an _ 1 a^n 
as ~ 2 

We denote by n"((5,5) the solution of this equation over the 
whole real axis -00 < S < 00, with the boundary condition 
that it vanishes at 6 = ±00. One can check immediately that 



(20) 



nV,5) = 



1 



-5-/(25) 



(21) 
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Fig. 2. — A trajectory that performs multiple up-crossings of the threshold 
at 5 = Sc- 



This probability distribution would bring us back to PS the- 
ory, and to its problems discussed in the Introduction. So, 
we need to eliminate the trajectories once they have reached 
the threshold. In B ond et al.l (1199 Ih this is implemented by 
imposing the boundary condition 



m,s)\ 



(5=5, 



= 0. 



(22) 



This seems very natural, but we stress that this boundary con- 
dition is still something that it is imposed by hand. The 
so lution of the FP equ ation with this boundary condition 
is (IChandrasekhadl943h 



m,s)= 



1 



-5-/(2S)_ -(25,-5) V(2S) 



(23) 



and gives the distribution function of excursion set theory. 
When studying halo merger trees it is important to consider 
also the distributi on for trajectories that start frora an a rbi- 
ti-ary value Sq ^0 jBond et al.lll99Tt iLacev & ColelfT993h . In 
this case, eq. (|23] ) is replaced by 

1 



n(^o;5;5): 



-(5-5o)-/(2S) _ ^-(25,-5o-5)7(2S) 



(24) 



This result is easily understood writing 26^ -d[)-S = 2(5^ - 
So)-(S-So), so eq. ( |24] | is obtained from eq. (l23T l performing 
the obvious replacement 6 ^ S-Sq, and also 6c 6c - 6o, 
which expresses the fact that, if we start from 6o, the random 
walk must cover a distance 6c - 6o to reach the threshold. 

In the excursion set theory the distribution 11(6, S) is defined 
only for 6 < 6c, so the fraction F{S) of trajectories that have 
crossed the threshold at "time" smaller or equal to S cannot 
be written, as in eq. (|5]l, as an integral from 6 = 6c to 6 = +oo. 
Rather, we use the fact that the integral of 11(6, S) from 6 = -oo 
to 6 = 6c gives the fraction of trajectories that at "time" S have 
never crossed the threshold, so 



F(S)=l- 



d6U(6,S). 



(25) 



Observing that 11(6, S) = Il"(6,S)-n%26c-6,S), we see that 
d6lf(6,S)+ / d6Il^(26c-6,S). (26) 

OO J-oo 



Since n^((5,5) is normalized to one, 

d6lf(6,S)= I d6lf(6,S). (27) 

OO J 5c 

For the last term in eq. (l26] t, we write 5' = 26c - and 

poo 

d6U"(26c-6,S)= / d6'n\6',S). (28) 

Thus, one obtains 

F(S) = 2j^ d6n\6,S) = siic(^-^ , (29) 

where i/ = 6c/a(M), and one recovers the factor of two that 
Press and Schechter were forced to introduce by hand. The 
probability of first crossing the threshold between "time" S 
and 5'+ii5' is given by !F(S)dS, with 



dF 
'dS 



\6''' 



dS ■ 



(30) 



This can be easily computed by making use of the fact that 11 
by definition satisfies the FP equation ( l20b . so 



dn 

'85 



6c 



(5=5, 



/2^53/2 



-5,-/(25) 



(31) 



Observe that, in (5 = 6c, Tl(6,S) and all its derivative of even 
order with respect to 6 vanish, while all its derivative of odd 
order with respect to 6 are twice as large as the value for the 
single gaussian (ISTT i. So, this first-crossing rate is twice as 
large as that computed with a single gaussian, which is an- 
other why of understanding how one gets the factor of two 
that the original form of PS theory misses. 

The halo mass function follows if one has a relation M = 
M(R) that gives the mass associated to the smoothing of 6 over 
a region of radius R. We discuss below the subtleties associ- 
ated to this relation, and its dependence on the filter function. 
Anyhow, once M(R) is given, we can consider f as a function 
of M rather than of S(R). Then \dF /dM\dM is the fraction 
of volume occupied by virialized objects with mass between 
M and M+dM. Since each one occupies a volume V =M/p, 
where p is the average density of the universe, the number of 
virialized object n(M) with mass between M and M + dM is 
given by 

dF 



dn p 

dM= — 

dM M 



dM 



dM, 



dn 
dM 



p dF 
M Is 



dS 



dM 



P 



,d\na 



M2 dlnM ' 



(32) 



(33) 



where we used S=a^. Therefore, in terms of the first-crossing 
rate T(S) = dF /dS, the function /((t) defined from eq. ([T]i is 
given by 

f(<j) = 2a^n<J^). (34) 

Using eq. ( [3TT l we get the halo mass function in PS theory 
(with the factor of two computed thanks to the excursion set 
theory). 



dn 
dM 



1/2 



-5;/(2(T-) 



p dlna 



M2 d\nM ' 



(35) 



This is the result given in eqs. ([T]i and (|2|i. 
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The crucial point is how to associate a mass M to the filter 
scale R. For the sharp filter in A:-space this is clear The mass 
associated to a spherical region of radius R and density p is 
M = (4/3)TrR^p. For the other filters there is no unambiguous 
definition. A possibility often used is the following. One first 
normalizes W so that its maximum value is one. Calling W' 
this new dimensionless filter, one can define the volume V 
associated to the filter as V = J d^xW', and M = pV. This 
procedure seems reasonable, but still it is somewhat arbitrary, 
since one might as well chose a different normalization for 
W'. For a gaussian filter, this gives V = {I-k^^I^R^. For a sharp 
A;-space filter, on top of this ambiguity, there is also the fact 
that such a volume is not even well defined. In fact, the k- 
space filter in coordinate space reads 



1 sinM — MCOSM 
H^sha.p-i(r,-R) = -^-^ , 

where u = r/R = kjr and r = |x-x'|, which gives 



SmM — MCOSM 



and 



y =47r / drp-W' =\2'itR^ 



/■OO 


smM 


/ du 


COSM 


JO 


M 



(36) 



(37) 



(38) 



The integral of sin m/m gives 7r/2, but the limit for A — > oo of 
the integral of cos u from m = to m = A does not exist. If one 
just sets it to zero, without much justification , one finds the 
result V = 6tt^R^ which is sometimes quoted dLacev & Cold 
Il993h . In any case, it is clear that it is difficult to give unam- 
biguous numerical predictions for the halo mass function with 
a filter different from the sharp ;c-space filter 

The standard practice in the literature is to use the PS mass 
function, which can be derived from excursion set theory but 
only if one works with a sharp fc-space filter, and at the same 
time to use M = {A/3>)nR^ p, which is only valid for a sharp 
It-space filter Of course this is not consistent and cannot be a 
good starting point for the inclusion of the non-Gaussianities, 
since one would attribute to primordial non-Gaussianities fea- 
tures in the mass function which are due, more trivially, to the 
filter function. 

In principle, one can determine the halo mass function with 
a tophat filter in coordinate space by performing a Monte 
Carlo (MC) realization of the trajectori es obtained from 
a Langevin equation with colored noise (iBond etal.llT99l 
[Robertson et al.ll2008t) . However, our final aim is to get some 
analytic understanding of the effect of non-Gaussianities on 
the halo mass function, and to this purpose we need a good 
analytic control of the effect of the filter, first of all in the 
gaussian case. 

3. PATH INTEGRAL APPROACH TO STOCHASTIC 
PROBLEMS 

3.1. General formalism 

We have seen that the computation of the halo mass func- 
tion can be reformulated in terms of a stochastic process. We 
now show how to compute the probability distribution of a 
variable evolving stochastically, in terms of its correlators. In 
this paper we limit ourselves to gaussian variables, while in 
paper III of this series we perform the generalization to arbi- 
trary non-Gaussian theories. 

Let us consider a variable 5{S) that evolves stochastically 
with "time" S, with zero mean ((5(5')) = 0. For a gaussian 



theory, the only non-vanishing connected correlator is then 
the two-point correlator {5{Si)5(S2})c, where the subscript c 
stands for connected. 

We consider an ensemble of trajectories all starting at 5o = 
from an initial position 5(0) = 5{), and we follow them for a 
time S. We discretize the interval [0,5] in steps A5 = e, so 
Sk = ke with k= and S„ = 5. A trajectory is defined 

by the collection of values {5\,. . . ,(5„}, such that 5(Sk) = 5k- 
There is no absorbing barrier, i.e. 5(S) is allowed to range 
freely from -oo to +00. The probability density in the space 
of trajectories is 

W^(5o; <5i , . . . , 5„-Sn) = {SoiMSx )-5i)... 5d(5(S„) - <5„)) , 

(39) 

where, to avoid confusion with the density contrast 5, we de- 
note the Dirac delta by So- In terms of W we define 

dSi... / dSn-lW{6o',5u . . . ,Sn-i,Sn;Sn) , 
00 J-oo 

(40) 

where S„ = ne. So, n£((So;(5;5) is the probability density of 
arriving at the "position" (5 in a "time" S, starting from So at 
time S{) = 0, through trajectories that never exceeded Sc- Ob- 
serve that the final point 6 ranges over -00 < S < 00. For later 
use, we find useful to write explicitly that 11 depends also on 
the temporal discretization step e. We are finally interested in 
its continuum limit, ne=o, and we will see in due course that 
taking the limit e ^ of He is non-trivial. 

The usefulness of is that it allows us to compute the 
first-crossing rate from first principles, without the need of 
postulating the existence of an absorbing barrier Simply, the 
quantity 

rSc 

dSIl,(6o;S;S) (41) 



gives the probability that at time S a trajectory always stayed 
in the region S < 5c, for all times smaller than S. The rate 
of change of this quantity is therefore equal to minus the rate 
at which trajectories cross for the first time the barrier, so the 
first-crossing rate is 

T(S) = - [ ^ d6dsn,(So;S;S) (42) 



(where ds = 9/95), just as in eq. ( |30] |. The halo mass function 
is then obtained from this first-crossing rate using eqs. ([T]i and 
(I34I 1. Observe that no reference to a hypothetical "absorbing 
barrier" is made in this formalism. We will discuss below 
how, and under what conditions, an effective absorbing barrier 
emerges from this microscopic approach. 

To express ll^(6o;5;S), in terms of the two-point correlator 
of the theory we use the integral representation of the Dirac 
delta 

Sd(x)= / ^e--^\ (43) 



2tt 



and we write eq. ( [39] l as 



dXl_ dX i A,<5; / „-i Et, A,5(S,)\ 

27r • ■ ■ 2^ ^ ^ ■ 

(44) 

Observe that the dependence on So here is hidden in the cor- 
relators of S, e.g. ((5^(5 = 0)) = Sq. It is convenient to set for 
simplicity Jo = in the intermediate computations, and it will 
be easy to restore it in the final results. For gaussian fluctua- 
tions. 
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as can be checked immediately by performing the Taylor ex- 
pansion of the exponential on the left-hand side, and using the 
fact that, for gaussian fluctuations, the generic correlator fac- 
torizes into sum of products of two-points correlators. This 
gives 



W(So;Su---,S„;S„)= / PA e'^- ^'^'^ ^'^^ , 
where 



dXi dX„ 
27r In 



(46) 
(47) 



and Si = S{Si). Then 



neiSo;S„;S„)-- 



d5i . . -dSn-i / VX 



exp< i^XiSi-^^XiXj{Si6j)r > 



(48) 



3.2. Gaussian fluctuations with sharp k-space fllter 

As we have seen in Sect. |2l the computation of the halo 
mass function in the excursion set formalism with sharp k- 
space filter can be reduced to a Langevin equation with a 
Dirac-delta noise. Therefore, we now study the case in which 
5 has gaussian statistics (so only the two-point connected 
function is non-vanishing) and obeys the Langevin equation 
(fTSl l with a noise r]{S) whose correlator is a Dirac delta, 
eq. (fT9] l. Using as initial condition So = 0, eq. (fTsT l integrates 
to 



dS'rj(S'), 







so the 2-point correlator is given by 



(49) 



{d(Sd6(Sj)),= / dS / dS'{r,(S)v(S')) (50) 
Jo Jo 
= min{Si,Sj) = emin(/,y) = eA/j . 

Denoting by W°™ the value of W when 5 is a gaussian variable 
and performs a markovian random walk with respect to the 
smoothing scale, i.e. satisfies eqs. ( fTSI l and (fT9] l, we get 



WS"'(6o;Su---,5„-S„) = 



dXi 
'2^ 



i=l 



(2^/^ (detl)'/2 ^^P E-^'-^^"')'/^' 



(51) 



Given that A,y = min(/, j), we can verify that A~' is as fol- 
lows: (A~'),-, = 2for! = l,...,n-l, (A"')„„ = 1, and (A"'),-,,+i = 
(A"' ),+! ,■ = -1, for i= 1, . . . ,n- 1, while all other matrix ele- 
ments are zero. Furthermore, det A = 1. As a result, we get 

1 



WS"(5o = 0;<5,,...,(5„;5„): 



(27re)«/2 



X exp < 



1 

"2^ 



(52) 



This expression takes a more familiar form using the identity 
2Si(Si-5M) = (dM-dif-(6li-Sf), together with EVA^Ii " 



Sf) = Recall also that eq. ( fSTT i assumed as initial con- 

dition 5o = 0. The result for So generic is simply obtained by 
replacing Si Si - So for all ; > 0. Then, for / > the terms 
(Si+i - Sif' are unaffected, while in the last term of the sum 
S\ ^ S\ — Sq. Thus, for So arbitrary, we get 



wn<5o;<5i,...,J„;5„) = 



exp ' 



{Si^i-SiY 



(53) 

Observe that W^"'(So\S\, . . . ,S„;Sn) dSi . . .dS„-\ is just the 
Wiener measure (see e.g. chapter 1 of lChaichian & DemichevI 
( 120011) ). From eq. (|53]l we see that 



W^'^iSo-^S, S„;S„) = ^,(S„-S„.i)W"''^(So; Si 



where 



^,(AS) = 



1 



exp 



(Ml 
2e 



,S„-\;Sn-i) , 
(54) 

(55) 



(27re)i/2 

Equation ( l54l l expresses the fact that the evolution determined 
by eqs. ( fTSl l and (fT9] l is a markovian process, i.e. the probabil- 
ity of jumping from the position (5„_i at time S„-i to the posi- 
tion S„ at time S„ depends only on the values of S„-S„-i = AS 
and on S„-S„-i = e, and not on the past history of the trajec- 
tory. Integrating eq. ( l54b over Si,..., Sn-i from -oo to Sc we 
get the important relation 

Uf"iSo;S„-S„)= I ' dSn.i<i,{S,-Sn-i)nf'{So\S„.i-Sn-i), 



(56) 

which generalizes the well-known Chapman-Kolmogorov 
equation to the case of finite Sc. 

4. DERIVATION OF EXCURSION SET FORMALISM FOR 
GAUSSIAN FLUCTUATIONS AND SHARP /f-SPACE 
FILTER 

We now want to derive, fro m our "microscop ic" approach, 
the excursion set formalism of Bond et alj ( 119911) . As we have 
seen in Section |2l the result of Bond et al. holds for gaussian 
fluctuations and sharp ^-space filter, working directly in the 
continuum limit, and reads 



nS,((5o;<5;5) = 



1 



-(5-5„)V(25) _ -(2&-&a-SfH2S) 



(57) 



We want to prove eq. ( fSTl ) using our definition of lie as a path 
integral over all trajectories that never exceed Sc. Beside being 
a starting point for the generalization to arbitrary filter func- 
tions and to non-Gaussian theories, the derivation of the ex- 
cursion set theory from first principles has an intrinsic interest. 
In fact, in Bond et alj (il991h this result is obtained by postu- 
lating that the distribution function obeys a FP equation with 
an "absorbing barrier" boundary condition Ii{So\ S\ S)\s=Sc = 0. 
While the fact that ne=o obeys a FP equation follows from 
eq. ( fTSl ). the absorbing barrier boundary condition is rather 
imposed by hand. As we already mentioned, in the liter- 
ature on stochastic processes it is known that, in the gen- 
eral case, the distribution functio n Ilf(So;S;S) does not satisfy 
any s imple boundary condition dvan Kampen & OppenheimI 
Il972t [Knessl 1986). It is therefore interesting to see how, 
in the gaussian case with sharp A:-filter, an absorbing barrier 
boundary condition effectively emerges from our microscopic 
approach. 
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We first show that in the continuum limit we recover 
eq. i5% . Then, we examine the finite-e corrections. As it 
turns out, these corrections have a non-trivial structure which 
is quite interesting in itself. Our main reason for discussing 
them in detail, however, is that they play a crucial role in the 
extension of our formalism to a generic filter function and to 
non-Gaussian fluctuations. 

4. 1 . The continuum limit 

To compute Ilf" by performing directly the integrals over 
(5i , . . . , S„-i in eq. ( |40l i. and then taking the limit e ^ is very 
difficult, since the integrals in eq. (|40] | run only up to 5c, and 
already the inner integral gives an error function whose argu- 
ment involves the next integration variable. 

A better strategy is to make use of eq. (|56] |. This relation 
expresses the fact that, for gaussian fluctuations and sharp k- 
space filter, the underlying stochastic process is markovian. 
We change notation, denoting 5„ = 5, 5„-5n-\ = A(5, and S„-i = 
S, so S„ = S+e. For fixed 6, we have dSn-\ = -d{AS), and 
eq. ( |56] | becomes^ 

poo 



Using eq. (I55l l we see that 



In the limit e 
order in e eq. 



Js-s, 

we have *e(A(5) 
gives 



(58) 

So(AS), so to zeroth 



/•oo 

nfo((5o;J;5)= / d{AS)SD{AS)nf!:^0Q-5-AS;S). (59) 

If 6- 6c < 0, the integral includes the support of the Dirac 
delta, and we just get the trivial identity that n^™Q((5o;(5;5) is 
equal to itself. However, if 6- Sc > 0, the right-hand side van- 
ishes and we get nJ™)((5o;(5;5') = 0. The same holds if 5 = Sc. 
In this case only one half of the support of is inside the in- 
tegration region, so we get n^™o('^o; S;S) = (1 /2)n^™)((5o; 5;S), 
which again implies n^™Q((5o; 6; S) = 0. Therefore we find that 

^l=o(So; S-S) = ifS> 5c . (60) 

This is not in contrast with the fact that nf™((5();(5;5) is the 
integral of the positive definite quantity For finite e, 

nf'"((5o;i5;'S) is indeed strictly positive but, when 5 > 6c, it 
vanishes in the limit e 0"^. 

Consider now eq. (ISST i when 5 < Sc- In this case the zeroth- 
order term gives a trivial identity. Pursuing the expansion 
to higher orders in e we have to take into account that in 
Uf^(5o;S;S + e) there is both an explicit dependence on e 
through the argument S+e, and a dependence implicit in the 
subscript e. We begin by expanding the left-hand side as 

dur(So;S;S) 



Iff (So; S;S+e) = nf™(5o; S;S) + e- 

d^uriSo;S;S) 



dS 



■2 dS^ ^''^ 
without expanding for the moment the dependence on the 
index e. On the right-hand side of eq. ( |58] ), we expand 
nr(So;S-AS;S) in powers of A^, 



d(AS) ■^AAS)Ilf^(So;S-AS; t) 



(62) 



s-s. 



=0 



(_l)"5«nr((5o;(5;5) 
dS^' 



d{AS){AS)"-^,(AS). 



5-5, 



' In this section we always assume that &q is strictly smaller than &c- The 
case (5o = &c is important when we study the non-markovian corrections, and 
will be examined in due course. 



d(AS)(AS)"^AAS)-- 



5-5,. 



(5,-5)/ 



dyy"e 



(63) 

If S is strictly smaller than Sc and Sc-S is finite (more pre- 
cisely, if it does not scale with ^/e) , the lower limit in the 
integration goes to -oo as e ^ 0^, and 



(5,-5)/V27 



dyy" 



dyy"e-^^ + O (e-('5"-«)V(2e) ] (54) 



1 +(-!)« 



(n-l)!! + C' 



2 2«/2 ^' V )- 

The residue, being exponentially small in e, is beyond any 
order in the expansion in powers of e, and we can neglect it, 
so 

1*00 

I d{A5){AS)"^AAS)'^f"'^(n-\)\\, (65) 
Js-s, 



if n even, and vanishes if n is odd. Thus, eq. ( I58l ) gives 

dIir{S,;S;S) S d^Iir{S,;S;S) 
il° (oo;o;o) + e — ^— ttt:^ 1" • 



dS 



-n, iSo,S,S) + - + -g + ... 

(66) 

From this structure it is clear that, when 6c-S is finite, the 
dependence on the index e in Hf™ can be expanded in integer 
powers of e, 

nnSo; S; S) = nf „(<5o; <5; S)+enf^^(So;S;S)+e^Uf^^(So; S;S)... 

(67) 

where H^", flf^, etc. are functions independent of e. We 
can now collect the terms with the same power of e in the 
expansion of eq. ( |66] |. To order e we find 



dIll%(So;S;S) I d^Ili:o(So-^d;S) 



= 0. 



(68) 



dS 2 dS^ 

Putting together this result with eq. (|60] l, we therefore end up 
with a FP equation with the boundary condition n^™y(5o;^ = 
Sc;S) = 0, and therefore we recover eq. ( |57] i. We have there- 
fore succeeded in deriving the excursion set formalism from 
our microscopic approach. Observe that the boundary con- 
dition nJ")((5o;^ = Sc\S) = emerges only when we take the 
continuum limit, and does not hold for finite e. 

4.2. Finite-e corrections 

In Section l53] we will find that the halo mass function gets 
contributions, that we will call "non-markovian", that depend 
on how nf™((5o; 5c;S) approaches zero when e ^ 0. It is there- 
fore of great importance for us to understand the finite-e cor- 
rections to the result obtained in the continuum limit. The 
issue is quite technical and we summarize here the main re- 
sults. Details are given in appendixlAl 

As long as Sc - (5 is finite and strictly positive, we have seen 
that the expansion (l67b applies, so the first correction to the 
continuum result is 0{e) and is given by ell^". Collecting the 

next-to-leading terms in eq. 



we find that n|^™ satisfies 
a FP equation with the second time derivative of nj"(, as a 
source term, 

dn^^^jSo; S; S) 1 d^f\^':^^{5o; 5;S) _\ d^If^^iSo; S; S) 
dS 2 dS^ ~4 dS^ 



(69) 
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In the above derivation, a crucial point was that we could ex- 
tend to -oo the lower integration limit in eq. ( l63T l. This is 
correct if we take the limit e — * 0""" with 6c - S fixed and pos- 
itive. The situation changes at (5 = S^, since in this case the 
lower limit of the integral is zero, rather than -oo. In this case 



fl'(AJ)(A(5)*,(A<5)= 



(70) 



(while the same integral computed from -oo to +oo obviously 
vanished), so we now have a term O (Ve) on the right-hand 
side of eq. i6% . Furthermore, 



^ 



so the expansion of eq. ( 1581 ) now gives 



e \i/2 dUr(SQ;5;S) 



(a 



as 



t=4 



This indicates that nr((5o;(5,.;5) is 0(e^/\ rather than 0(e). 
However, eq. ( |72] i is not a good starting point for a quantitative 
evaluation of Ilf^{do;Sc',S) since, as we show in appendix lAl 
the expansion in derivatives becomes singular in 6 = S^, and all 
terms denoted by the dots in eq. (|72] | finally give contributions 
of the same order in e. A better procedure is the following. 
First, observe that the correction is determined by the lower 
limit of the integral, (Sc - S)/ V2e. The transition from the 
behavior 0{e) valid for Sc - S fixed and positive, to the behav- 
ior ©(e'/^) valid at S = Sc takes place in a "boundary layer", 
consisting of the region where Sc-S is positive and ©(e'/^), 
and the lower limit of the integral is 0(1). This is a situation 
that often appears in stochastic processes near a boundary, or 
in fluid dynamics, and can be treated by a standard technique 
(see e.g. Knessl (1986), where a very similar situation is dis- 
cussed in terms of the means first passage time, rather than in 
terms of the distribution function Ilf™). Namely, we introduce 
a "stretched variable" 77 (not to be confused, of course, with 
the noise ry(r) of eq. (fTsT l) 



77: 



Sc-S 



(73) 



which even as e ^ 0""" is at most of order one inside the bound- 
ary layer, and we write nf™(i5o; S; S) in the form 



Uf^{So;S;S) = C,(So;S;S)u(7^), 



(74) 



where Ce(So;S;S) is a smooth function, while the fast variation 
inside the boundary layer is contained in uirj). By definition, 
we choose u{ri) such that lim,,_^oo "(J?) = 1, so is just the 
solution for Ilf" valid when Sc -S is finite and positive, i.e. 
Cf is given by eq. (|67] |. Writing S = Sc-rj^/le (and setting for 
notational simplicity So = 0) we have 



CASo = 0-S;S)=- 



1 



(75) 



exp 



Sc-riV2e 



-exp 



1 / 



plus corrections 0{e). Since Q by definition is smooth every- 
where, we can use eq. ( iTST i also inside the boundary layer In 



this case t] is at most 0(1), and we can expand the exponen- 
tials in eq. ( 1751 ) in powers of y/e. In the limit e — > 0, 

Plugging this result in eq. ( |74] i and sending S ^ S' we find 



where 



nf"(5o; Sc;S) = ^/i ^^^-^^^Ps) + ^^^^^ ^ 



-= Urn r}u{rj) . 

/TT I)— ►0 



(71) In appendix lAl we show that 7=1/ ^/tt, so 



(77) 
(78) 

(79) 



(72) Similarly, for Sn < Sc, 



UriSc;S„;S) = V~e4= ^-^e<'^-'"^' 1'-^'^ 



+ 0(e). (80) 



0F 53/2 

Observe that at the numerator of eqs. ( |79l ) and dSOb always 
enters the absolute value of the difference of the first two ar- 
guments of nf™, i.e. Sc-S{) in eq. (l79b and Sc-S„ in eq. (l80b . 
as it is also obvious from the fact that Ilf™ is definite pos- 
itive. Equations (|79] l and dSOl l will be important when we 
compute the non-markovian corrections, in Section 15.31 To 
conclude this section, it is interesting to discuss the behavior 
of nf "((Jo; S; S) for S larger than Sc, with S-Sc finite (and, as 
always in this section, Sq < Sc). In this case the lower integra- 
tion limit in eq. ( 1631 ) goes to +00 as e ^ 0"^ and 



nr(So;S;S): 



.exp{-iSc-Sf/i2e)}. (81) 



This function is zero to all orders in a Taylor expansion around 
e = 0+. 

5. EXTENSION OF EXCURSION SET THEORY TO GENERIC 

FILTER 

We next consider the computation of the distribution func- 
tion He, still restricting for the moment to gaussian fluctua- 
tions, but using a generic filter function. In this case the nat- 
ural time variable is the variance S computed with the cho- 
sen filter function, so in the following S denotes the vari- 
ance computed with the filter function that one is considering. 
Again we discretize it in equally spaced steps, Sk = ke, with 
S„ = ne = S, and a trajectory is defined by the collection of 
values {Si ,S„}, such that S{Sk) = Sk- 

The distribution function for gaussian fluctuations and arbi- 
trary filter function is given by eq. (|48] |. As we saw in the pre- 
vious section, in the markovian case lie satisfies a local differ- 
ential equation, namely the Fokker-Planck equation. It is in- 
structive to understand that, for a generic filter, it is no longer 
possible to write a local diffusion equation for Ile(SQ;S„;S„). 
This will immediately make it clear that the problem is now 
significantly more complex. Indeed, by taking the derivative 
with respect to Sn of both sides of eq. ( |48] |. we get 



lle(oo;o„,i„) = } / 



dSn 



(82) 



X / dSi...dS„-idkdiW(SQ;Si,...,S„;S„), 
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where dk = d/dSk, and we used the fact that, acting on 
exp{/^"^j XjSi}, dk gives iXk- Therefore, separating the term 
with k = I = n from the rest, and observing that {S{Sk)S{Si))c 
depends on S„ only if at least one of the two indices ^ or / is 
equal to n, we get 



d 1 92 

oSn 2 dof, 



(83) 



+i] ^^^^9„ dSi . ..d5n-i dkW{5a; Su---, 5nA). 



k=l 



If the upper limit of the integrals were +00, rather than Sc, the 
term proportional to dkW with k <n would give zero, since it 
is a total derivative with respect to one of the integration vari- 
ables dSi ,.. . d5„-\ , and W vanishes exponentially when any 
of its arguments 5k goes to ±00. Thus, one would remain with 
a Fokker-Planck equation. However, when the upper limit 5c 
is finite, the terms proportional to dkW with k < n give in 
general non-vanishing boundary term. Actually, for a sharp 
A:-space filter, we found that {5k5„)c = mm{Sk,S„) = Sk, which 
is independent of S„ for k < n. Therefore d{5k5„)c/dS„ = 0, 
and the term in the second line of eq. (|83] | vanishes. This is 
another way of showing that, in the continuum limit, for sharp 
A:-space filter the probability distribution satisfies a FP equa- 
tion, as we already found in Section|4T|^ 

For a generic form of the two-point correlator, the term in 
the second line of eq. ( l83T l is non-vanishing, and in general it 
is very complicated. Furthermore, in the continuum limit the 
sum over k in eq. ( |83] | becomes an integral over an intermedi- 
ate time variable Sk, so this term is non-local with respect to 
"time" S. Thus, we can no longer determine Il^{5{);5„;S„) by 
solving a local differential equation, as we did in the marko- 
vian case. Once again, this shows that the common procedure 
of using the distribution function computed with the fc-space 
filter, and substituting in it the relation between mass and 
smoothing radius of the tophat filter in coordinate space, is not 
justified. What we need is to formulate the problem in such a 
way that it becomes possible to treat the non-markovian terms 
as perturbations, which is not at all evident from eq. ( l83T l. 

In this section we develop such a perturbative scheme. We 
illustrate the computation of Il^{5o;5„;S„) using a tophat fil- 
ter in coordinate space, which is finally the most interesting 
case since we can associate to it a well defined mass, but the 
technique th at w e develop can be used more generally. 

In Section lsTl we study the two-point correlator with tophat 
filter in coordinate space and we show that it can be split into 
two parts, which we call markovian and non-markovian, re- 
spectively. In Section 15.21 we compute the contribution of 
the markovian term to the halo mass function, while in Sec- 
tion 15.31 we develop the formalism for computing perturba- 
tively the contribution of the non-markovian term. 

5.1. The two-point correlator with tophat filter in 
coordinate space 

We first study the correlator {5{Ri)5{R2)) with a tophat filter 
in coordinate space. We use eq. (O with x = 0. The two-point 
correlator of the non-smoothed density contrast is given in 
eq. (fTTl i. We write the power spectrum after recombination 
as P(k)T^(k), where P{k) is the primordial power spectrum 

^ Note however that this only holds in the continuum limit, as is implicit 
in the fact that we are taking the derivative with respect to S„, which means 
that we are considering 5„ has a continuous variable. 
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Fig. 3. — The functions S(R) computed for a tophat filter in coordinate 
space (blue soHd curve) and S(R) computed for sharp i:-space filter (violet, 
dashed), against R, on a log-log scale. 



and T(k) is the transfer function, so for the smoothed density 
contrast we get 

1 f°° 

{5(R,)5(R2)) = / dkk^P(k)T\k)W(k,RiW*{k,R2). 
27r^ Jo 

(84) 

When Ri = R2 = R, this reduces to S{R). We consider a 
primordial spectrum P(k) = Ak"', processed into the post- 
recombina t ion sp ectrum by the transfer function T(k) as in 
ISugivamal (fT99l . in a concordance ACDM model with a 
power spectrum normalization = 0.8 and h = 0.7, flM = 
1 -17a = 0.28, Qb = 0.046 and «s = 0.96, consistent with the 
WMAP 5-years data release. 

We first study S{R). We compute the integral in eq. ( fTsT i 
numerically, for different values of R, both with the sharp k- 
space filter (|7]) with kj = l/R, and with the tophat filter in 
coordinate x-space, whose Fourier transform is 

sm(kR)-kRcos(kR) 

Wsharp-.(^, R) = 3 . (85) 



{kpy 



For both filters, the constant A in P{k) is fixed so that S 
when R = (8//1) Mpc. The result is shown in Fig. [3] 

We consider next the correlator ( l84l l with the tophat filter in 
coordinate space. We compute the integral in eq. ( |84| ) numer- 
ically, holding R2 fixed and varying Ri. The result is shown 
in Fig.|4] The solid line is the function S{Ri), already shown 
in Fig.m The dashed fine is {6(Ri)6(R2)) with R2 = 1 Mpc/h, 
as a function of Ri, while the dotted line is {5(R])5(R2)) with 
R2 = 5Mpc//!, again as a function of We see that, as long 
asRi < R2, the two-point correlator is approximately constant 
and equal to S{R2), while for 7?i > R2 the correlator is approx- 
imately equal to S{Ri). In other words, 

{5(Ri)5(R2)) ~ mm(S(Ri),S(R2)) ■ (86) 

In Fig. |5] we compare {d(Ri)d(R2)) (blue solid line) and 
min{S(R]),S{R2)) (violet dashed line). This result suggests 
to define a function C{Ri,R2) from 

{6(Ri)6(R2)) = nim(S(Ri),S(R2)) + C(Ri,R2). (87) 

As we see from Fig. [3j the function S{R) can be inverted to 
give R = R(S), so Ri = R(Si) and R2 = R(S2). We define 

A(Si,S2) = C(R(,Si),R(,S2)), (88) 
and we can write 



{5i5j)=min(Si,Sj) + A(Si,Sj). 



(89) 
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Fig. 5. — The coiTelator {5(R[)5(R2)) (blue solid line), compared to 
min(S(Ri),S(R2)) (violet, dashed). In both cases R2 = IMpc/h, and we plot 
the functions against R[, measured in Mpc//i. 



If one simply neglects A(5',,5';), i.e. one makes the approx- 
imation (l86b . the problems become formally identical to the 
one that we have solved in Section H] Therefore, we end up 
with the standard excursion set theory result given in eq. ( l23T l, 
and therefore with the PS mass function, in which S is sim- 
ply the variance computed with the filter of our choice, in this 
case the tophat filter in coordinate space. The corrections to 
this result are due to A(Si,Sj), so it is useful first of all to 
understand better the form of this function. 

By definition A(Si,Sj) is symmetric, A{Si,Sj) = A{Sj,Si), 
so it is sufficient to study it for 5, < Sj. We also use the no- 
tation Ajj = A(Si,Sj). Since, by definition, (Sj) = Si, we see 
from eq. ( [89l ) that A(Si,Sj) vanishes when S, = Sj. Further- 
more, at Si = 0, (5, = S{) is the same constant for all trajectories, 
so {SiSj)c = So{Sj)c = 0, and therefore A{Si,Sj) vanishes when 
Si = 0. 

In Fig.|6]we plot A(Si,Sj) for Sj fixed, as a function of Si, 
with <Si < Sj, for our reference ACDM model (solid line). 
The dashed line in Fig.|6]is the approximation 



MSi,Sj)'. 



,Si(Sj-Si) 



(90) 



with K ~ 0.45 (a more accurate value will be given below). 
We see that eq. ( |90l ) provides an excellent analytical approxi- 
mation to A(5;,5/).^ 

' Varying Sj we find that eq. 190) becomes exact (within our numerical 



Fig. 6. — The function A(S,,5j) for tophat filter in coordinate space (solid 
line), with Sj ~ 1.42 (corresponding to R(Sj) = 5 Mpc/ft), plotted against 



5;, in the range < < Sj 
(dashed line). 



and the function KSi{Sj — Si)/Sj with k ~ 0.45 



For Sj fixed and 5, 0, the correction A(Si,Sj) is linear 
in Si, so more generally we can define K{Sj) from K(Sj) = 
lims;^o A{Si,Sj)/Si, or equivalently. 



k(R)= lim 



{S(R')S(R)) 



1. 



(91) 



S'^oo {SHR')) 

In the ACDM model that we are using, our numerical results 
display a very weak linear dependence of k on R. Taking for 
instance the data in the range R E [l,60]Mpc//i, the result of 
the numerical evaluation of eq. (IgTT l is very well fitted by 

k(7?)~ 0.4592-0.0031/;, (92) 

where R is measured in Mpc//!.'° 

5.2. Markovian term 
we get 



Inserting eq 

ne((5o;(5„;5„) = 



into eq 
dSi 



. dS„- 



vx 



(93) 



X exp A/i 



\,-,5,- + ^(^i, S j)] XiXj 

i=i i.j=\ 

As we see from Fig. |5] eq. ( l86l l gives a reasonable approxi- 
mation to the exact correlator. This suggests to treat A,y as a 
perturbation, so we now expand in A,j. The zeroth-order term 
is simply ni™((5o;5„;5„), whose continuum limit is given in 
eq. ( |57] |. The corresponding first-crossing rate is 



1 



2^ 53/2 ' 



dS 



(94) 



accuracy) for small Sj, while for large Sj the function Si(Sj—Si)/Sj must be 
replaced by a a less symmetric function, whose maximum is at a value of 
5, slightly larger than Sj/2. The qualitative shape of the function remains 
however the same. For completeness, we have also considered a gaussian 
filter In this case we find that, in a first approximation, the function ^ij 
is still given by eq. )90t (although the actual form of A,; is slightly more 
skewed compared to an inverse parabola), with a value of k ~ 0.35. 

The value of k depends in principle on the cosmological model used, but 
this dependence is quite weak. For comparison, using a ACDM cosmological 
model with h = 0.7, = 1 - f^A = 0.3, erg = 0.93, UbIi^ = 0.022 and n,, = 
1, consistent with the WMAP 1st year data release, gives k(R) ~ 0.4562 — 
0.0040R. 
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so the markovian term can be obtained by taking the excursion 
set result dSTb . which was computed with the sharp fc-space 
filter, and replacing the variance computed with the sharp k- 
space filter with the variance computed with the filter of inter- 
est. This is the procedure that is normally used in the litera- 
ture. From our vantage point, we now see that the corrections 
to this procedure are given by the non-markovian contribu- 
tions, to which we now turn. 

5.3. Non-markovian corrections 

We now discuss the non-markovian corrections, to first or- 
der, using the analytical approximation ( |90l l for A,y. From 
eq. (|93] ). expanding to first order in A,j and using A,e' ^» = 
-idie' ^1 , where 9,- = d/dSi, the first-order correction to 
is 

nfi (5o;<5„;5„)= / dSi . . .d6„.i -J^^iAdj 

X y"l?Aexp|/^A,,5,-i^min(5,,5,)A,A,| (95) 

1 " r^' 

= 9 E ^'V / d5,... dS„.i didjW^'^iSo; <5i , . . . , S„;S„) . 

We rewrite the term Aijdidj separating explicitly the deriva- 
tive dn = d/dS„ from the derivatives 9,- with / < «, so (using 

- n - n—l 1 

-Y,^udidj=-A„„d^„ + J2AM + -J2^ijdidj. (96) 

i,j=l i=l i.j=\ 

Since A,y = when / = j, the above equation simplifies to 

^ n n-\ 

- ^ijdidj = E ^'"^'■^« + E ' ^97) 



where 



1=1 /<j 

n-2 n-l 



E-EE 

i<7 /=1 j=i+l 



(98) 



When inserted into eq. ( |95] l the term X^Li' ^mdid„ brings a 
factor that, in the continuum limit, produces an integral 
over an intermediate time 5,. Because of this dependence on 
the past history, we call this the "memory term". Similarly, 
the term X];<j ^U^'^j gi^^^, in the continuum limit, a double 
integral over intermediate times S, and Sj, and we call it the 
"memory-of-memory" term. Thus, 

jj A 1 _ jjmem _|_ jjmem-mem 



where 

nr'"(5o;<5„;5„) 



(99) 
(100) 



= E A,>,a„ / dSi... d5„-i diW^^'iSo; Su...,S„-S„), 



and 



n mem-mem/ c . c . c \ 
J (do;d„;S„) 



(101) 



= E A,7 / t/^i . . . dS„-ididj WS'"((5o; 5 1 , . . . , 5„) . 



If we expand to quadratic and higher orders in A,y, we get 
terms with a higher and higher number of summations (or, in 
the continuum limit, of integrations) over intermediate time 
variables. 

To compute the memory term we integrate 9,- by parts, 
dSi... dS„-i 9,WS'"((5o; Si,...,Sn; S„) (102) 

-oo 

d5\ . ..d5i...d5„-iW{5Q;5\,.. .,6/ = Sc,.. . ,S„-uS„;S„) , 

where the notation dSj means that we must omit dSj from the 
list of integration variables. We next observe that, because of 
the property ( l54l i. W^™ satisfies 

WsySo; , . . . , <5,-i , (5,, ,...,S„;S„) (103) 

= W^'^iSo; Su..., Si-i , S,; SdW^'^iS,; ,...,S„;S„- Si) , 

so 

dSi . . . dSi-i / dSj+i . . . dS„-i 

-oo J —oo 

X s,^..., Si-i , Sc-, SdW^'^iS,; , . . . , 5„ ; S„ - Sd 

= nf^(6o-6,.;SdIlf^(Sc;S„-S„-Sd , (104) 
and we get 

«-i 

nr"(^o;<5„;5„) = EA«,9„ [nr(So;Sc;Sdnr(6c;S„;S„-Si)] ■ 

(105) 

In the continuum limit we write 



E 



1 



dSi , 



(106) 



and, using eqs. ( |79] l and dSOl l. we find 

n»7((5o = 0; 4; 5„) =-d„ f" dSi A(5,-, S„) ~ ^"^ 

TT 



X exp 



s]'^{S„-Sd^l^ 
5] {5c-5„f 



(107) 



25,- 2(5„-5,). 

We now insert the form ( |90l l for A,y. The integral can be 
computed exactly using the identities 



s„ 



dSi ■ 



Si 



■ 5p(5„-5,)^V2 I 25,- 2(5„-5,) 
/—I 1 r (a + fc)^ 



(108) 



and 



S„ o2 

s^'\s„-Siy^^ 

,1/2 



25,- 2(5„-5,) 



where Erfc is the complementary error function. ' ' This gives 
Hf:^{So = 0;S„;S„)=Kdn 



SciSc-S„) f 26c-S„ 
Erfc 



Sn 



/2,Sn 



(110) 



' ' To derive these results we take one derivative of the left-hand side of 
eq. 1108) with respect to a^. The resulting integral can be performed using 
eq. )A5t . and we then integrate back with respect to a^. Similarly, eq. 1109) 
is obtained taking twice the derivative with respect to a^. 
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as functions of o. 




FIG. 8.— The functions nf™(5o = 0; 5;S = l)+n™™™™(<5o = 0; 5;S = 1) 
(dashed), compared to n^™Q((5o = 0;(5;5 = 1) (solid line), as functions of 5. 

For the memory-of-memory term, proceeding as for the mem- 
ory term, we get 

n™'"-"™((5o;(5„;5„) (111) 
= Ai.nr (<5o; 5c\ 5,)nr(5,; 6,; Sj - 5,)nr (5,; 6„-S„ - Sj) . 

To compute this quantity we also need Ilf"{Sc',Sc;S), with 
both the first and the second arguments equal to Sc. As we 
discuss in appendix|A] the result is 

nn^,;4;5)=^=^. (112) 



Actually, eq. (II 12l i is exact, and not just valid to 0(e). Using 
eqs. ( |79] l and ( II 12b we get 

:(5,(4-<5„) 



TTinem- 


mem 








dSi 












dSj 







ttx/Itt 



(113) 



-exp 



(Sc-S„f 



It is convenient to use the identity 

(Sc-S„) ^ f (Sc-S„f \ f (Sc-S„f 

(S„ - Sj) ^''P 1 2(S„ - Sj) J "^""^X 2(S„ -Sj)^' ^ 



to write n™™"""^™ as a total derivative with respect to 6„ . The 
inner integral can now be computed rewriting it in terms of 
the variable z = iSc-S„)^/[2(S„-Sj)], and gives 



n'T"'"™('5o = o-dn;S„) = -^=^ 9„ 

\/2ttS„ 

X ^e-''/^^''^Erfc ({Sc-S„) 



-{S-S„r/(2S„) 







Si 



Si 



.(115) 



We have not been able to compute analytically this last in- 
tegral, but the fact that n™™"™™ is a total derivative with 
respect to S„ will allow us to compute analytically the first- 
crossing rate, see below. First, it is interesting to plot the 
functions Hf™ and n^o™"™™. We show them in Fig.Hl set- 
ting for definiteness S„= 1 . Observe that these two functions 
are separately non-zero in 5 = Sc. However, 



n'i7(<5o = 0;(5,,5„) = -^Erfc 



and n™5"-™'=™ 



(116) 



{6,, = 0; 5c,S„) = -n™7(Jo = 0; Sc,S„), so we find 



that the total distribution function Il^^)(5o',S;S„) still satisfies 
the absorbing barrier boundary condition Ilg=o{6[);Sc;Sn) = 0, 
even when we include the markovian corrections to first order. 
In Fig. [8] we compare 11 
term ( |57t 



'"r + n™7"'°™ to the zeroth-order 



5.4. The halo mass function 

We can now compute the first crossing rate using eq. ( [30l ). 
Since both 11™™ and n™™"™'=™ ha ve be en expressed as a 
derivative with respect to d„ in eqs. dl 10b and dl 151) . the in- 
tegral over d5„ is performed trivially, and we get 

dS. 



■mem-mem 



(S)-- 



fl'5„n™7(<5o = 0;^„;5) = 0, 
(So = 0;6„;S) 



(117) 



-OS I 



mem-mem / 

£=0 



d_ 
"dS 



■oo 
kSc 



1 



2'kS Jq Si 



kSc d 



V2^dS 
kSc 



51/2 [ '2S 



(118) 



2S^n r'25 



where r(0,z) is the incomplete Gamma function. Putting to- 
gether eqs. ( |94] |. dl 17b and dl 18b we find the first-crossing rate 
to first order in the non-markovian corrections, 

ns) = ^ j^g""-^^''^ + r 1 0, ^ 1 . (119) 



f2i sy^^ 

The halo mass function in this approximation is therefore 



f{a) = il-K)[- 

TT 



1/2 



20-2 



2tt <J 



(120) 

where, in the relevant range of values of R, k is given by 
eqs. d9Tl ) and d92l ), and is a slowly decreasing function of R. 
For instance, at 7? = 5 Mpc, k ~ 0.45, at = 10 Mpc, k ~ 0.43, 
and at = 20 Mpc, k ~ 0.40., For large values of 5c /2a^ 



(121) 
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Fig. 9. — The function f{u) against u. The gray dotted Hne is the value in 
PS theory. The blue solid line is the fit to the numerical W-body simulation 
of Warren et al. 1 2006.) . The thick black dashed line is our result (120). 



Thus the incomplete Gamma function gives the same expo- 
nential factor as PS theory but with a smaller prefactor, so for 
large halo masses it is subleading, and eq. ( I120l i approaches 
(1 - k) times the PS prediction. 

In Fig.|9]we plot the function fiiy), where i^ = 6c/(t, compar- 
ing the prediction o f PS theory given in eg. (|2|i, the fit to t he A^ - 
body simulation of IWarren et alJ (l2006l) . and our result (I120ll. 
This figure can be compared to Fig. 4 of [Robertson et al.l 
(EPOS), see in particular their bottom-left panel, where the 
authors show the prediction of PS theory, the result of their 
A^-body simulation, and the computation of f{iy) with tophat 
filter in coordinate space, performed with a Monte Carlo real- 
ization of the trajectories obtained from a Langevin equation 
with colored noise. We have used the same scale and color 
code as their Fig. 4, to make the comparison easier One sees 
that our analytical result for f{iy) agrees very well with their 
Monte Carl o result (the function th at we call f( i') is denoted 
as i^fiv) in [Robertson et al.l (l2008l) ). From eq. (1120b . we see 
that in the end our expansion parameter is just k, so evalu- 
ating the non-markovian corrections to second order we will 
get corrections of order k^. For k given by eq. (|92] | these are 
expected to be of order 20%, which is the level of agreement 
between our analytical result and the Monte Carlo computa- 
tion. This provides a non-trivial check of the correctness of 
our formalism. 

A second consistency check is obtained by recalling that 
the fraction of volume occupied by virialized objects is given 
by eq. ( l25l l. In hierarchical power spectra, all the mass of the 
universe must finally end up in virialized objects, so we must 
have F{S) = 1 when 5c/(t ^ 0. Formally, the limit 5c/(t 
can be obtained sending (5c — > for fixed a, so we require that 

Um / d5IV{5,S) = Q). (122) 

As we recalled below eq. (|5]l, the original PS theory fails this 
test, giving that only one half of the total mass of the universe 
collapses. Inourcase n = ns'" + n™™ + n"™-"™. Since Hs™ 
is the same as in the standard excursion set result, it already 
satisfies eq. (1122b . so we must find that, in the limit 5c ^ 0, the 
integral of n™"" + 11™'='""™'=" from -oo to 5c vanishes. Using 



eq. (1107b we see that 

(123) 

for all values of Sc- For the memory-of-memory term we find 

£,,na~(W) = ^|jr(o,|). (124, 

Since, forz — + 0, r(0,z) -Inz, we have 

so eq. (1122b is indeed satisfied. An equivalent derivation starts 
from the observation that, in terms of the function /(cr), the 
normalization condition reads 

da 

/ — /(a)=l. (126) 

Jo cr 

Substituting /(cr) from eq. ( |120b into eq. ( |126b and using 

da(^"^ ^e-''^l^^^'-^=l (127) 

and 

we see that the dependence on n cancels and eq. ( |126b is 
satisfied. The term proportional to the incomplete Gamma 
function therefore ensure that the mass function is properly 
normalized, when the amplitude of the term proportional to 
exp{-(5^/(2CT^)} is reduced by a factor 1 - k. 

A number of comments are now in order. First, our findings 
confirms the known result (Bond et al. 1991; Robertson et al] 
l2008h that the corrections obtained by taking properly into 
account the tophat filter in coordinate space do not alleviate 
the discrepancy of PS theory with the A^-body simulations. 
We see in fact from Fig.|9]that the effect of the non-markovian 
corrections is to give a halo mass function that, in the relevant 
mass range, is everywhere smaller than the PS mass function, 
which results in an improvement in the low-mass range but in 
a worse agreement in the high-mass range. This indicates that 
some crucial physical ingredient is still missing in the model. 
This is not surprising at all since, as we already stated, the 
formation of dark matter haloes is a complex phenomenon. 
Incorporating some of the complexeties within the excursion 
set theory will be the subject of paper 11. 

On the positive side, we conclude that we have developed a 
powerful analytical formalism that allows us to compute con- 
sistently the halo mass function when non-markovian effects 
are present. In this paper we have applied it to the corrections 
generated by the tophat filter function in coordinate space. 
However, the same formalism allows us to compute pertur- 
batively the effect of the non-Gaussianities on the halo mass 
function. This direction will be developed in paper 111. 

Before leaving this topic we observe that, in the perturba- 
tive computation performed in this section, all terms turned 
out to be finite in the continuum limit. The fact that the to- 
tal result is finite is obvious for physical reasons. However, 
the fact that all the terms that enters in the computation are 
separately finite happens to be a happy accident, related to 
the form (|90] l of A(Si,Sj), and in particular to the property 
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A(5', ,5',) = 0. However, not all perturbations that we will con- 
sider share this property, in particular when we consider the 
non-Gaussianities. Furthermore, even with the above form of 
A{Si,Sj), if we work to second order in perturbation theory 
we find that divergences appear It is therefore important to 
understand in some detail how in the general case these di- 
vergences cancel among different terms, giving a finite result. 
This issue is quite technical, and is discussed in detail in Ap- 



pendix B. 
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APPENDIX 

A. FlNlTE-e CORRECTIONS 

In this appendix we derive the results for n(:((5o;(5c;'5) and Il^(Sc;Sc',S) mentioned in Section l4!2l These results are needed in 
Section l531 when we compute perturbatively the non-markovian corrections. We will also show that, for dsHf^iS^ = 0; 6; S), the 
limit e ^ 0"^ does not commute with the limit 6 6~. This result will be important in appendix|B] when we study the cancellation 
of divergences that can appear in intermediate steps of the computation. In eqs. ( iTTl i and (iTSl l we found that nf"((5o; (Jc;^) = 
^/e ^iSc/S^^^)e~^'-^''^^^ + 0(e), where 7 = (2/^/t:) lim^^()?7M(77). One possible route to the evaluation of 7 could be to plug eq. ( |74] | 
into eq. ( fSST l and evaluate both sides at 6 = Sc- To lowest order in e one can replace 5 + 6 on the left-hand side simply by S, and 
one obtains an integral equation for the unknown fu nction ^(77). This integral equation has the form of a Wiener-Hopf equation, 
for which various techniques have been developed (" NobldflQSS !). However, we have found a simpler way to get directly 7, as 
follows. We consider the derivative of Ilf™ with respect to 6c (which, when we use the notation ni™(5();(5;5), does not appear 
explicitly in the list of variable on which Ilf" depends, but of course enters as upper integration limit in eq. (|40l)). This gives 

^n<:((5o;4;5„) = / d6i...d6i...d6n-iW(do;6i,...,6i = 6c, ■■ .,6^,8,1), (Al) 
d6c Z^J-°o 

where the notation d6i means that we must omit d6i from the list of integration variables. We next use eqs. ( 1103b and ( 1104b and, 
in the continuum limit, we obtain the identity 

^n,=o(<So;5„;5„)= / "t/^,- lim-nf'"(5o;<5,;5,)nf'"(5,;,5„;5„-5,). (A2) 
d6c Jo e^oe ' 

The left-hand side of this identity can be evaluated explicitly using eq. (ISTb and, setting for simplicity 60 = 0, is 

±IIM^ = 0;6„;S„)= (-V' ,«)V(2^,.) . (A3) 

d6c \ttJ 53/2 

The right-hand side of eq. (IA2b can be evaluated using eq. ( |77] ) together with 

n,(4;4;5) = n,(<5„;5,;5) = V^7^|^e-(^'-^">'A2S) + 0(g)^ (A4) 

which can be checked from eqs. (|40] i and (|53] l by performing a reshuffling of the dummy integration variables. We see that the 
limit e — > in eq. ( IA2I ) is finite thanks to the factors in Ilf"(6o',6c;Si) and in Ilf"(6c;6„;S„- Si). The integral over Si can be 
performed using the identity 

where a > Q,b > 0.^^ In this way we find that the dependence on 6c and S on the two sides of eq. ( IA2b is the same, as it should, 
and we fix 7 = 1 / ^/tt. 

In appendix|B] when we study the cancellation of divergences, we will also need dsUf^, evaluated in 6 = 6c- Of course, if we 
first take the limit e ^ 0"^, and then we take 6 6~, we simply get the derivative of the function ne=o(i5o;'5,5') given in eq. (|57] |. 
evaluated in 6c, 

/2\l/2 r 

£rn £rn^55nr(<5o = 0;(5,5)= a,-nf:„(^o;5,5)|,^, =-f-j ^e"*'/^'''- (A6) 

However, we will actually need the result when the limits are evaluated in the opposite order, i.e. lime^o+ lims^s- dsIlf"{6Q; 6,S). 
We will now show that these two limits do not commute. From eq. (fTsT l. for small 77, u{ri) is proportional to jiT-j]) = 1 /{2rj). 
More generally, for small 77, we write 

1 1 

u{r])=— + UQ + Uy7^ + 0{rf). (A7) 
Zrj 

We have not been able to find this identity in standard tables of integrals, but we have verified it numerically, with very high accuracy, in a wide range of 
value s of a and b. We can also turn the argument around and say that, sinc e we know that n|"'(<5o;<5r;S) has the functional form (77) and we know that the identity 
)A2t holds, it follows that the integral on the left-hand side of eq. )A5t must be given by the expression on the right-hand side, times an unknown numerical 
constant. The latter can be computed evaluating the term ~ 1/a of the integral in the limit a 0^. This is easily done analytically, since in this case the factors 
(S„ — Si) inside the integrand can be simply replaced by 5„, and fixes the factor \/2tt on the right-hand side of eq. jA5). 
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Plugging this expansion, together with the expansion in powers of 77 of eq. (fTST l. into eq. (|74] | we find that, for 77 ^ (i.e. for 
S ^ 6~ at fixed e), retaining only the terms up to Oiy/e) 

Uf-iSo = 0; S, S) = nr((5o = 0;5,,S) + V~e^ iu^v + ^d + ■ ■ ■)^^~'''^^'' ■ (A8) 

Using ds = {dri/dS)d/ drj and drj/dS = -1/ \/2e, this gives 

/2\i/2 r 

hm lima5n?"(<5o = 0;<5,5) = -Mo - ^e"*'/*'", (A9) 

which differs by a factor mq from eq. ( IA6I ). It is also interesting to observe, from eq. ( |A8t , that also 9^nf"/977^, evaluated in 
77 = 0, is proportional to ^/e. Since 

92nr_i52nr 

overall S^nf™/^^^, evaluated in (5 = 5c at finite e, is proportional to l/y^. Therefore, in eq. ( |72] | the first correction included 
in the dots, which is proportional eid^Iif^ / d6^)s=s^., is of the same order as the term ^/t{^Iif^ / d6)s=5,, and similarly for the 
higher-order terms. This is the reason why we could not use eq. (|72| i to fix the value of the coefficient 7. 

Finally, in the perturbative computation we also need nf™((5c;(5c;5), with both arguments equal to 5c- The result is given in 
eq. ( 11121 ). To derive it, we first observe from eq. ( |79] l that, when 5q = 5c, the term 0(\/e) vanishes, so the first non-vanishing term 
will be Oi f). Invariance under space translations requires that nf"((5o; i5c; S) can depend on 5o and 5c only through the combination 
5c - (5o, so when 5q = 5c it becomes a function of S only. We can perform dimensional analysis assigning to 5 some (unspecified) 
dimension I and to S dimensions f-. In this case, from eq. (fT9] l we see that 77 has dimensions l/£, and ^ ~ i/l'^ = l/£, so eq. ( fTsT i 
is dimensionally correct. In these units A ~ 1 /£, since X5 is dimensionless, and we see from eq. (1481 ) that has dimensions 1 /£. 
Using dimensional analysis in this form we conclude that the term 0{e) in ni™((5e;(5c;5) is necessarily proportional to tjS^I'^. 
Since this fixes completely the dependence on 5, writing S = €nv^t have also fixed completely the dependence on n, i.e. to 0{e) 
we must have 

nn<5.;<5.;5) = c^ = -^, (All) 

with c independent of n. The coefficient c can then be fixed computing explicitly the integral in eq. ( |40] l when n = 2, i.e. when 
there is just one integration variable. This can be done analytically and shows that c = 1 / \/27r. The computation for n = 2 actually 
shows that eq. ( lAl 11 1 is exact, i.e. it receives no correction of higher order in e. Even for n = 3 the integral in eq. ( lAl 11 1 can be 
performed analytically when 5q = 5c, and again we find that eq. (lAl lb is exact. We have checked this result numerically for n up 
to 7 and we fi nd tha t the numerical result agrees with eq. (lAl lb within the 10 digit precision of the numerical integration, so it is 
clear that eq. ( lAl lb is actually exact, and not just the result at 0(e). (In any case, to perform our perturbative computation, we 
only need nr(5c;5c;S) to 0{e).) 

B. DIVERGENCES AND THE FINITE PART PRESCRIPTION 

In this appendix we first of all reconsider the perturbative computation of Section l531 for a generic fu nctio n A,j (still symmetric 
in (j,7)). This will reveal some complexities that were not apparent in the computation of Section 1531 and that will be very 
important when computing the non-Gaussianities. If A,j does not vanish when / = j, we rewrite eq. ( |93] l as 

n,(<5o;<5„;5„)= / ' d5,...d5„-x / PAexp <( -> ' A; V exp <( / > 'A/^,--) '[min(5,-,5,)]A,A/ ) , (Bl) 



5„-i jvXexpl ^^Aijdidi >exp<| A/^,- i^[min(5,-,5^)]A,Aj > , 



where, as usual, we used the fact that, acting on exp{iA,(5,}, 9,- gives /A,-. Since A„„ is now in general non-vanishing, in the sum 
I the term A„„9^ contributes. Furthermore, now 

, «— 1 , n-l 

- ^ A,,9,9, = ^ A,jd,dj + A/,9f . (B2) 



2 ■ 



The operator exp{(l /2)A„„d^,} can be carried out of the integral over d6i, . . . ,d5n-i, while the other terms A,j will be expanded 
perturbatively. Thus, eqs. (|99b-( fTon i are replaced by 



where 



and 



jjAI _g(l/2)A„„9- l^-j-jmem _|_ jjmem-memj (B3) 

nr™(<5o; 5„-S„) = A,„a„ / d5i... d5„-i a,WS™(5o; <5i , . . . , 5„; 5„) , (B4) 

.■ 1 J -oc 

5]A„9,.9, + if]A,a; 



Ilf~'"(5o;5n;S„)= I d5,...d5„.i 



WS"'{5o;di,...,5„;S„). (B5) 
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The memory term is the same as in Section 15.31 so it is still finite. The memory-of-memory term, however, presents a new 
difficulty. Using eq. ( 1103b we get 

n~™(^o; s„-s„) = ^u^nSo; 4; 5,)nf-(5,; 6,-Sj - SdnriS,.; s„-s„- Sj) 

A 

[^r(So;Si;Si)Ur(Sr,S„;S,-Sd] (B6) 

n-l 

= Ao-nr (,So; s.-Sdnf^iSc-, Sc;Sj - 5,)nr s„-s„ - Sj) + Y A,[9,nr(,So; Sr, Sd]5,=5M"'(Sc\ Sn;S„ - Sd . 

i<j 1=1 

We now discover that the continuum limit of the memory-of-memory term is non-trivial, since it is made of two terms that are 
separately divergent. Consider first the second term in eq. ( IB6I ). which is the one coming from A„9?. We have found in Section l4~2l 
that nf™((5o;(5f.;5,) is proportional to while [dsUf"i6c;S;S)]s=s, has a finite limit for e — > 0, see eq. iA9i . Therefore, using 
eq. ( 11061 ). we find that the last term in eq. (IB6I 1 diverges as 1/ y/e. A similar problem appears in the term coming from djdj with 
/ y j. Using eqs. ( |79] l and ( II 12l i we find that the first term in eq. (IB6b is proportional to 

where 5, = ie, Sj = je. In the continuum limit, unless A,y vanishes for / = j, this quantity diverges as 1 / -/i, because of the behavior 
(Sj-Si)~^^'^ when Sj ^ Sf. In Section l4~2l these problem did not show up because A,, = 0, so the divergence coming from A„9? 
disappears. Furthermore, when Si S j, A^ vanished as Sj-S,, thereby ensuring the convergence of the sum (or, in the continuum 
limit, of the integral over Sj) in eq. ( |B7| i. 

In order to understand how the cancellation mechanism works when A,y does not vanish for St = Sj, we examine the memory- 
of-memory term when A{Si,Sj) is a constant, that we set equal to unity. The reason is that, in this case, we can compute it in an 
alternative way, which shows that the result is finite. The trick is to compute the second derivative of Ilf" with respect to Sc. The 
first derivative was computed in eq. (lAlb . and the result can be rewritten as 

d r^' 

—nf^(So;6„;S„) = V / dSi... d6„-i diW . (B8) 

When we take one more derivative of eq. (lAlb with respect to 5c, we find two kinds of terms. First, there are the terms where we 
take one more derivatives with respect to the upper limit of the integration with respect to a variable d5j with j ^ i. Furthermore, 
we must take the derivative of W{6{)\ (5i , . . . , (5,- = (5^, . . . , (5„-i , 5„;S„) with respect to 5c. Therefore 

^nr(^o;5«;5„) = 2^ r dS^. ..TSi . ..T5j . ..d5„-i W(5o;5u. ..,5i = 5c, ...,5j = 5c,.. .,5„;S„) 

^"e j^j J-oc 



(B9) 



that is, 



Thus, when A,j = 1, 



/■ ' ---- d 

+ / d5\...d5i...d5„-\—W{5o;5i,...,5i = 5c,...,5„;Sn) 
jz^J-oo dSc 

d5i...d6„-ididjW + Y / d5i...d6„-idfw , 

,<,j °° i=i 

—^Uf^(5o;d„;S„) = J2 d5i . . .dd„.ididjW . (BIO) 

i j=l -00 

U~-'-(So;5„;Sn)=^-^Iir{do;5„;S„). (Bll) 



In particular, in the continuum limit. 



nSr"™(<5o = 0;5„;S„) = i^nf!'o(<5o = 0;5„;S„) = (^^ ' 



^_ (15c- 5„f 



Sft 



1 ^-(25,-5„)V(2S„)_ 



First of all this result shows that, when A,y = 1, n™™"™™ stays indeed finite in the continuum limit. Second, it gives its explicit 
expression, which can then be compared with a computation based on eq. ( IB6b . To perform the comparison, we first compute the 
second term in eq. ( |B6b . when A,j = 1, i.e. 

n-l 

h = J2[diUf^i5o = Q;5r,Si)h=snf"(Sc;5„;S„-Sd. (B13) 
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Observe that in this expression we must first compute the derivative in (5, = Sc (since this came from the integration by parts of 
df) and only after we take the limit e 0""". The result is therefore given by eq. ( IA9l l. Using also eqs. (IA4b and ( |79] |. we get 



1 U0V2 ^ 1 r iSc-s„y- ^ ^^^^^ 

Ve I- 5p(5„-5,)3/2 I 2Si 2(5„-5,), 



Because of the exponential factor, the argument of the sum goes to zero very fast as 5, 0^ and as 5, 5„, and therefore we can 
use eq. ( 11061 ). so 

1 U0V2 , , 1 r iSe-S„f \ 

The integral can be performed using eq. ( IA5b . and we get 

/, =-4 ^(25,-5„) (B16) 

Therefore /] diverges as 1 / It is important to observe that there is no finite part in /j . In the continuum limit the corrections to 
eqs. (|79] l, (|A9l l and (1106b are all 0{e) compared to the leading terms that we used, so they produce terms that are overall 0{^/e) 
in eq. ( IB16I ). and therefore vanish in the continuum limit. 
We next consider the other term in eq. ( IB6I ). i.e. 



n—2 n—\ 

l2 = J2 nrC'^o = 0; 6,; 5,) ^ 0^(4; S,; Sj - 5,)nr(<5,; (5„; 5„ - Sj) 
i=i j=i+i 

n-2 «— 1 



. _L..,,,-.,.g ^.-'^'"".g -p{-tg} ■ 

Now the passage from the sums to integrals is more delicate. One might be tempted to write 

= / dSj. (B18) 

j=i+\ ■'^i 

However, eq. ( IB 18b is only correct when the sum and the integral are finite for e ^ 0. Here this is not the case, since 

diverges at the lower integration limit Sj = S,, and indeed our aim is to extract this divergent term, plus the finite terms. A better 
guess would be that, since the sum starts from j = i+l, the corresponding integral should start from Sj = 5, + e, so 

Still, this cannot be completely correct. To realize this observe that, since the integral is dominated by Sj = 5, + e, the divergent 
part can be extracted replacing Sj = Si everywhere except in the factor (Sj - Si)~^/^ so, if we used this prescription, we would 
conclude that 



h = TTT exp < > / dS: TT7 + finite parts 



exp < > + finite parts . (B2 1) 



V^(5„-5,)3/2 2{S„-Sd 

However, if the prescription ( IB20I ) where correct in general, we should get the same result if we separate the term j = i+1 from 
the sum, and we let the remaining integral start from Sj = S, + 2e, so we should get the same result if we write 

, ? 1 1 r (5c-s„f \ , ,„ 1 r iSc-s„)^ 

h = —j= TTT exp < r / " '^J TF) TF) ^^'^ 



V?(5„-5,)3/2 2(5„-5,)J 4+2e ■'(5,-5,)3/2(5„-5,)V2 2(Sn-Si) 

exp < > + finite parts . (B22) 



^/i (5„-5,)3/2 2(5„-5,) 

We see that the two procedures both agree on the fact that the singularity is proportional to 1 / \/e, but give different values for 
the coefficient, so eq. (IB20b cannot correct in general. Observe also that the finite parts are not affected by this ambiguity, which 
amounts to a rescaling of e. 
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Since, of course, the strength of the singularity is in principle fixed (although difficult to compute analytically) as long as we 
write h as a sum, we can always choose a value a such that, as far as the 1 /\/e singularity and the finite terms are concerned, we 
have the equality 



«-i 

j=i+i 



1 



{Sj-SifiKs„-Sj)^i^ 



exp 



i5c-5nf 



dS 



1 



Si-van 



exp 



'2{S„-Sj) 



(B23) 



and the two expressions only differ by terms that vanish as e ^ 0. In fact, a can be fixed requiring that the coefficient of 1 / 
is the same on the two sides of eq. (IB23I ), and it does not affect the terms 0{€^) since it it just a rescaling of e. Actually, in our 
problem, an even better way to pass from the sum to the integral is to write 



n-l 

j=i+\ 



1 



{Sj-Sd^HSn-Sj)^l^ 



exp 



2{S„-Sj) 



dSj 



1 



{Sj-Si)yHs„-Sj)y^ 



exp 



ae 



2(Sj-Si) 2(S„-Sj) 



(B24) 



In other words, rather than setting the integrand to zero for Sj < Sj + ae, we cut it off exponentially using the factor exp{-ae/(Sj- 
Si)}. Again this produces a 1/ ^/ae singularity, as we will check in a moment, and a can be chosen so that this singularity has the 
same strength as that on the left-hand side of eq. ( IB24b . However, since a is just a rescaling of e, a does not affect the finite part. 
The advantage of using eq. ( IB24b is that the integral can now be performed analytically using eq. ( IA5b . so we get 



n-l 

./■='■+ 



1 



'^(Sj-sdyHSn-Sj)y^ 



exp 



'2(S„-Sj) 



1 



1 



1 



exp 



2(5„-5,) 



(B25) 



We have also checked this result numerically. The sum on the left-hand side can be computed very easily numerically, say for 
n up to 10'*, and we find that the right-hand side reproduces it perfectly, for all values of (5„, Si and if we choose a ~ 0.92. 
Expanding the dependence on at in the exponential and omitting the terms that vanish in the limit e ^ 0, we find 

1 1 1 ____ f (<5c-<5„)' 



h 



1 



Sn Si J 5c 5„ 



1 

(S„-Si)y^ 



exp 



(Sc 

'2(Sn-Si) 



(B26) 



which explicitly displays the 1 /y/e singularity and the finite, a-independent, part. 

To compute I2 we must still plug this expression into eq. ( IB 17b and carry out the sum over /. The latter sum presents no 
difficulty since its argument converges well both at Si = and at Si = S„, so we can just replace the sum by an integral using 
eq. ( |106b . It is actually convenient to leave 73 in the form eq. ( IB25I ), so the integral over 5, can again be performed using eq. iA5h 
and we finally get 



/2 = 



1/2 



1 

73/2' 

On 



-(2Sc-5„f/(2S„) 



2dc-S„ 



(2Sc-5„f 



S„ 



Putting together this result and eq. ( IB16b we finally find 



n mem— mem 
e=0 



(So;S„;S„) = 



-Mo 



V2] i25c 



On 



-(2(5,-(5„)7(2S„) _ 



1/2 



(2Sc-6„f 

Sn 



On 



-(2S,-5„)-/(2S„) 



(B27) 



(B28) 



However, in this case we already know the exact result for 11™™ ™'^™, which is given by eq. (IB 12b . Comparing these two results 
we learn the fo llowi ng. First, we know from eq. (IB12b that the result is finite and there is no l/\/e term. In the computation 
lead ing t o eq. ( IB28b we rather find two separately divergent contribution, so they must cancel. This is fully consistent with 
eq. ( IB28b . since these divergent terms have exactly the same dependence on e, 5^, 5n and Sn- We also see that, in this second way 
of performing the computation, the cancellation depends on the numerical values of quantities, such as mq, that are determined by 
the solution in the boundary layer, and which therefore are difficult to compute, as well as on the constant a that we determined 
numerically. The finite part is instead completely fixed, and it is not affected by the solution in the boundary layer, nor by the 
constant a, and correctly reproduces eq. ( IB12b . 

From this explicit example we can now extract a general rule of computation. Whenever A(Si,Sj) is a regular function, such 
as that given in eq. ( l90l l. the memory-of-memory term and analogous quantities which are finite when A{Si,Sj) = 1, will still be 
finite. The explicit computation with the formalism developed in Section [53] can generate terms that are separately divergent 
when e — > 0"*". However, since the total result is finite, these divergences must cancel among them. When we find integrals that 
diverge in the limit in which two integration variables become equal (such as the limit Sj Si above) we can just regularize 
them as in eq. (IB 24b . We call this technique "the a-regularization". We then discard the divergence and we extract the finite part, 
which is independent of a. We will indicate by the symbol J-V this procedure of taking the finite part. In this notation, the result 
of the above computations can be summarized by 



;=1 



n—2 n- 1 

^^EE 

,= 1 j=M 



dSi . . .d6„-\didjW ■■ 



1/2 r 



{25c- 5nf 
S„ 



On 



-(2(5,-<5„)-/(2S„) 



(B29) 
(B30) 
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As an application of the above formalism, we have studied what happens choosing a different expansion point when computing 
the halo mass function with a tophat filter in coordinate space. Observe in fact that, since A{Si,Sj) is symmetric under exchange 
of Si with Sj, eq. ( |90l l. which is valid for Si < Sj, can be rewritten more generally as 



MSi,Sj): 



min(Si,Sj)- 



[mm(Si,Sj)] 
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(B31) 



max{Si,Sj) 

Thus, the two-point correlator can be written as 

{Si6j) = (l + K)mm(Si,Sj) + A(Si,Sj). (B32) 

where, for 5, < Sj, A{Si,Sj) = -nSf/Sj. We can therefore use (1 + K)mm{Si,Sj) as the unperturbed two-point function, and treat 
A,y as the perturbation. The zeroth order term can again be computed exactly, since it just amounts to a rescaling of S, 5— > + 
At first sight this seems to give a modified exponential in the distribution function, since factors such as exp{-(5^ / (25')} in eq. dSTT i 
becomes exp{-S^ / [2{1 + k)S]} . However, now A„„ = -kS„ is non-zero, and we should not forget the factor exp{(l /2)A„„9,^} in 
eq. ( |B3t . The effect of this term can be computed exactly using the identity 

exp I \(b-a)d;'j -L ^ ^ ^-.V(2.) ^ (B33) 

which is valid for a > and b > 0. To prove it, we write 

dXsr^ ]_ f tZ^y (i\\2niXx-a/2)aX' _ f°°^-(l/2)(b-a)\'i\x-(l/2)a'\^ 



2^ «' V 2 ; 27T 

^ iXx-dmbX' _ J_ -.r/(2fe) /o^4N 

2^ Vb ■ ^ ' 

(Observe that for b <Q the final integral over dX does not converge, so this identity only holds if b > 0). In this way, we find 
that the action of exp{(l/2)A,„,c),^} on exp{-S^ / [2(1 + k)S]} gives back the "unperturbed" exponential factor exp{-(5^/(25)}, 
so the zeroth-order term of this expansion is finally the same as eq. dSTl i. The computation of the non-markovian corrections 
requires the finite part prescription, since now A(Si,Sj) does not vanish for 5, = Sj. The integrals over dSi and dSj are more 
difficult to compute, but for 6^ / (25) 3> 1 their exponential dependence is easily computed and, after taking into again the operator 
exp{(l /2)A„„(9^} in eq. (IB3b . we find that the exponential dependence of the corrections is the same that we obtained in eq. (II 19b . 
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